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ABSTRACT 

We analyzed the statistics of subhalo abundance of galaxy-sized and giant-galaxy-sized halos 
formed in a high-resolution cosmological simulation of a 46.5Mpc cube with the uniform mass 
resolution of lO^Af0. We analyzed all halos with mass more than 1.5 x lO^^M© formed in this 
simulation box. The total number of halos was 125. We found that the subhalo abundance, 
measured by the number of subhalos with maximum rotation velocity larger than 10% of that 
of the parent halo, show large halo-to-halo variations. The results of recent ultra-high-resolution 
runs fall within the variation of our samples. We found that the concentration parameter and the 
radius at the moment of the maximum expansion shows fairly tight correlation with the subhalo 
abundance. This correlation suggests that the variation of the subhalo abundance is at least 
partly due to the difference in the formation history. Halos formed earlier have smaller number 
of subhalos at present. 

Subject headings: cosmology: theory — galaxies: dwarf — methods: n-body simulations 



1. Introduction 

The cold dark matter (CDM) model is widely 
accepted as the standard theory of the forma- 
tion and evolution of the universe. According to 
this CDM model, the structure formation in the 
universe proceeds hierarchically, from small-scale 
structures to larger-scale ones. 

The CDM model successfully reproduced large- 
scale structures (e.g., Davis et "aLlll985 : Springel et al 
20051 ). However, it has been claimed that there 
are serious discrepancies between the structure 
predicted by the CDM model and that observed, 
in the scale of galaxy-sized objects and below. 
Moore et al] (|l999a[) and iKlvpin et all (|l999[ ) 



calculated the evolution of galaxy-scale dark mater 
halos using high-resolution cosmological iV-body 
simulations. They found that dark matter ha- 
los of the mass comparable to the Local Group 
contained far too many subhalos compared to 
the number of known dwarf galaxies in the Local 



^National Astronomical Observatory, Mitaka, Tokyo 
181-8588, Japan ; ishiyama@cfca.jp, makino@cfca.jp 

^Department of General System Studies, College of Arts 
and Sciences, University of Tokyo, Tokyo 153-8902, Japan 

^K&F Computing Research Co., Chofu, Tokyo 182- 
0026, Japan, fukushigOkfcr.jp 



Group. In the case of dark matter halos of the size 
of typical clusters of galaxies, the theoretical pre- 
diction and observation agreed pretty well. The 
theoretical prediction for galaxy-sized and cluster- 
sized halos are similar, and observations of clusters 
of galaxies and that of Local Group are very dif- 
ferent. 

This discrepancy, now called "missing dwarf 
problem" , was studied by ma ny researchers . 
Follow-up simulation studies (e.g., Dien iand et al 



i2004; Reed ct al. 2005; Kase et al. 20 03 gave es- 



sentiall y the same results as tho se of iMoore et al 



(|l999al ) and iKlvpin et al.l (|l999t ). Accordingly, it 
is now regarded as one of the most serious prob- 
lems of the CDM model. 

A number of solutions for this problem have 
been proposed. For example, we can reduce the 
number of subhalos in small scales by chang- 
ing the nature of dark matter in small scales. 
Models in this direction are warm dark matter 
(jKamionkowski fc Lidd"i3l2000l) and self-interacting 
dark matter ( Spergel &: Steinhardt 2000l ). An- 
other direction is to explain the discrepancy as 
the difference between the number of nonobserved 
dark matter subhalos and that of observed satel- 
lite galaxies. Maybe not all subhalos host dwarf 
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galaxies, or the measured velocity dispersion of 
dwarf galaxies might be much smaller than the 
rotation velocity of the parent subhalos. Models 
in this direction includes the s uppression of star 
form ation by early reionization (jSusa fc Umemura 
20041 ) and self-regul ation of st ar formation in sma ll 
halos (jStoehr et all 120021 : iKravtsov et all 12004 . 
However, none of these explanations are widely 
accepted as the clear-cut solution for the missing 
dwarf problem. 

In almost all previous studies of substructures 
in dark-matter halos, the so-called re-simulation 
method has been used. In this method, one sim- 
ulation is done in the following two steps. In the 
first step, a relatively large volume (for galaxy- 
sized halos typically a 50-lOOMpc cube) is simu- 
lated with low-mass resolution, and the candidate 
regions for high-resolution simulations are identi- 
fied. In the second step, the selected regions are 
simulated with higher mass resolution. In practi- 
cally all recent high-resoluti on simulations, this re- 
simulation method is u s ed (iDiemand et al.l 2007 , 
20081 : IStadel et al.ll2008l : fSpringel et al.ll200j^ . 



This re-simulation method allows us to resolve 
the structure of a selected halo with a very high 
resolution. On the other hand, it is not clear what 
kind of selection bias is introduced by the use of 
this method. Typically, the authors claimed that 
the halo which looks similar to the halo of the Lo- 
cal Group was selected to compare with the Local 
Group. However, in most cases the selection cri- 
teria were not described in detail, and thus it is 
difficult to see if any selection bias affected the 
results or not. 

An obvious way to avoid the selection bias, 
which might be introduced by the use of the re- 
simulation method, is simply not to use it and 
do the high-resolution simulation of the enti re ini- 
tial simulation box. In llshivama et all l|200al . here- 
after Paper I), we simulated a relatively small re- 
gion (the size of the simulation box is 21.4 Mpc) 
with high-mass resolution (mass of particles = 
3 X IO^Mq). We identified 21 dark halos with 
maximum rotation velocity larger than 200kms^^ 
and investigated the environmental effect on the 
subhalo abundance. The main conclusion of Pa- 
per I is that the abundance of subhalos has large 
variation. The most subhalo-rich halos contain 
the comparable num ber of subhalos as reported in 
Moore et al.l (jl999al ) , while least subhalo-rich ones 



contain around 1/10 subhalos. They also found 
that only the galaxy-sized halos (mass less than 
3 X IO^^Mq) show large variation in subhalo abun- 
dance, and more massive ones show much smaller 
variation. The subhalo abundance shows corre- 
lation with the concentration parameter c (larger 
c, smaller abundance) and measures of the local 
density such as the smoothed local density and 
distance to the nearest massive halo. 

The accuracy of the result of Paper I is, how- 
ever, limited by both the mass resolution (parti- 
cle mass 3 x 10^ Mq) and small number statistics. 
The total number of halos we analyzed was 21, and 
only 10 had the mass less than 3 x IO^'^^Mq. In this 
paper, we report the result of much larger, and 
thus much more accurate, simulation than that in 
Paper I. We simulated a box of the size 46.48 Mpc 
with 1600^^ particles of mass of 1 x 10^ Mq. We 
also reduced the softening parameter from 1.8kpc 
to 700 pc. Thus, compared to the calculation in 
Paper I, we used three times better mass resolu- 
tion, two times better spatial resolution, and 10 
times larger simulation volume. We identified all 
halos with the maximum rotation velocity larger 
than 160kms~^ and the virial mass M in the range 
of 1.5x1012Mq <M< IxIO^^Mq. We found 125 
such halos. In all halos, the number of particles in 
the virial radius is more than 1.5 x 10^, which is 
large enough to give reliable results for the statis- 
tics of subhalos with rotation velocity more tha n 
1/10 of that of the parent halo (jKase et al.ll2007l ). 

In Section 2, we describe the initial model and 
the simulation method. In Section 3, we present 
the result. Section 4 presents summary and dis- 



2. Initial model and Numerical Method 

We adopted a LCDM {^o = 0.3, Aq = 0.7, 
h — 0.7, as — 0.8) cosmological model for all of 
our calculations. In the largest simulation (sim- 
ulation I), we used a cube of the comoving size 
of 46.48Mpc and periodic boundary condition as 
the simulation box. The number of particles is 
1600'^ and the mass of one particle is 10^ Mq. To 
generate initial particl e distributions, we u sed the 
MPGRAFIC package (|Prunet et al.ll2008l) . which 
i s a parallelized var iation of the GRAFIC package 
(jBertschingerll2001l ). The initial redshift was 65. 

For the time integration, we used the GreeM 
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code (T.Ishiyama et al. in preparation), which is 
a massively parallel Tre ePM code based on the 
parall el TreePM code of lYoshikawa fc Fukushige 
()2005l) . It uses the recursive multisection algo- 
rithm (i Makino 2004 ) modified to achieve better 
load balancing for domain decomposition and the 
Phantom GRAPE sp ecial force calculation code 
( Nitadori et al.ll2006[ ) modified for the force with 
cutoff (Nitadori & Yoshikawa in preparation). The 
number of the grid points for PM calculation was 
400'^. The opening parameter of the tree part was 
= 0.30 up to z = 10, and was 6 = 0.50 from 
z = 10 to z = 0. 

We integrated the system using a leapfrog in- 
tegrator with shared and adaptive time steps. The 
step size was determined as min(2.0-\/e/|3i|, 2.0e/\vi 
(minimum of these two values for all particles). 
The total number of time steps was 15,844. The 
(plummer) gravitational softening e was constant 
in the comoving coordinate up to z = 10, and was 
constant (700 pc) in the physical coordinate from 
z = 10 to z = 0. 

We used 2048 CPU cores on the Cray-XT4 at 
Center for Computational Astrophysics, CfCA, of 
National Astronomical Observatory of Japan. The 
calculation time per step was about 60 s. The full 
run (15,844 time steps) was completed in 280 hr 
(wallclock time). 

In order to investigate the effect of the spatial 
resolution, we performed two additional simula- 
tions (simulations SI and S2). We used the same 
initial condition for these two runs and changed 
the softening parameter. In table[Tl we summarize 
the number of particles A'', the box size L(Mpc), 
the softening £(pc), the number of steps Ngtep, and 
the redshift Zg after which the softening was con- 
stant in the physical coordinate, for the three sim- 
ulations. 

The method to find halos and subhalos is the 
same as that used in Paper I. Our method is based 



Table 1: Run parameter 



Parameter 


Sim I 


SI 


S2 


N 


1600^ 


800^ 


800^ 


L(Mpc) 


46.48 


23.24 


23.24 


e(pc) 


700 


700 


350 


A'stcp 


15844 


10278 


20513 




10 


10 


21 



on the idea of finding all local potential minima. 
The identification of the main halos is done in the 
following three steps. In the first step, we sam- 
pled particles randomly according to the sampling 
rate /isamp, and calculated their potential. In this 
stage, all sampled particles are candidates for the 
centers of halos. Next, we select the particle with 
the smallest (most negative) potential and list it as 
the center of a halo. We then remove nmin neigh- 
bor particles of this center of halo particles from 
the original list of all particles, and select the po- 
tential minimum particle from the remaining par- 
ticles. At this time, we again look at nmin neigh- 
bor particles from the list of originally selected 
particles, and if one of the neighbors has smaller 
potential, we do not add this particle to the list 
of halos. However, we remove nmin neighbors no 
matter if the particle is added to the list of halos 
or not. We repeat this procedure until there are no 
remaining particles. Next, we calculated the halo 
virial radius i?v in which the sphe rical overdensit; 
is 178f2o'^ times the critical value ( Eke et al]|l99 
, the halo virial mass M and the maximum rota- 
tion velocity Vp. For identifying the main halos, 
we set i?samp = 0.01 and nmin = 1000 (this value 
corresponds to 1.0 x IO^^Mq). 

The method of finding the subhalos is almost 
the same as that for parent halos. The difference 
is that we do not reduce the particles. We detected 
all local potential minima inside the halo virial ra- 
dius i?v as subhalos. For identifying the subha- 
los, we set nmin = 100 (this value corresponds to 
1.0 X 10*Mq). 

3. Result 

We analyzed all halos with the maximum rota- 
tion velocity larger than 160 kms^^ and the virial 
mass 1.5 x IO^^^q < M < 1 x lO^^^^ ^he ro- 
tation velocity is defined as the maximum value of 
Vcirir) = y/GM{r)/r . We found 68 galaxy-sized 
halos (1.5 X 10^2^/0 < M < 3 x IQ^'^Mq), and 
57 giant-galaxy-sized halos (3 x lO^^M© < M < 
1 X lO^^M©). 

3.1. Effect of the spatial resolution 

Let us first discuss the reliability of our result. 
Since the current calculation used relatively large 
softening, the difference between the results of SI 
and S2 runs indicates the numerical error of the 
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Fig. 1. — Number of subhalos of SI plotted against 
those of S2. Here, Nx,i,Sj is the number of sub- 
halos for halo i with normalized rotation veloc- 
ity larger than x, in simulation Sj. Open cir- 
cles, squares and triangles show iV>o.i5, -^>o.io 
and iV>o.o5, respectively. 

Table 2: The difference of the subhalo abundance 



X 


0.15 


0.10 


0.05 


< N^.t^si > 

< Nx,i,S2 > 

< AN^/Nx > 


14.4 
16.5 
0.175 
0.246 


45.0 
54.8 
0.233 
0.258 


183.9 
220.4 
0.190 
0.200 



main calculation. We analyzed 21 galaxy-sized 
halos and 13 giant-galaxy-sized halos in simula- 
tions SI and S2. In these two runs, we changed 
only the softening and other parameters are all 
kept the same. The reason why we test the ef- 
fect of the softening here is that we u se the soft 



ening significantly larger than that in iKase et al 
(I2OO7I). thoug h com parable to what is used in 



Springel et al. (|2008l) for comparable mass reso- 
lution (simulation Aq-A-5). 

It is known that subhalos contain less than 200 
particles is not reliable (Kase et al. 2007'). The 
subhalo abundances iV>o.i5, Af>o.io satisfy this 
condition. However, -/V>o.o5 does not. Here, N^x 
is defined as the number of subhalos with the max- 



imum rotation velocity larger than x times that of 
their parent halos. 

Figure [T] shows the subhalo abundance of indi- 
vidual halos in run SI plotted against those for 
the same halos in run S2. Open circles, squares, 
and triangles show iV>o.i5, iV>o.io, and A^>o.o55 re- 
spectively. We can see that halos in run S2 with 
smaller softening have more subhalos than their 
counterparts in run SI. The difference, however, 
is small. 

Table [2] shows the average relative difference be- 
tween the number of subhalos in run SI and that 
in run S2, and their standard deviation. They are 
calculated as 



(1) 



■halo i \ jVi.i.Sl 



where rihaio is the number of halos and Nx,i.Sj is 
the number of subhalos for halo i with normalized 
rotation velocity larger than x, in simulation Sj. 
We can see that the deviation between the results 
of runs SI and S2 is smallest for x — 0.05. This is 
partly due to small number statistics. As shown 
in tabled! the number of subhalos is much smaller 
for higher rotation velocity. 

The value of A'o.os is, however, not reliable be- 
cause the number of particles in subhalos with this 
rotation velocity is too small. Therefore, we use 
the value iVo.i as the measure of the subhalo abun- 
dance. It does contain the error due to large soft- 
ening, but as we have seen in this section that is 
small. So we ignore this error in the analysis given 
in this paper. 

3.2. Subhalo abundance 

Figure [2] shows the cumulative numbers of 
subhalos as a function of their maximum ro- 
tation velocities Vc normalized by those of the 
parent halos at z = 0. We can see that the 
variation of the subhalo abundance is signifi- 
cantly larger for galaxy-sized halos (the top panel) 
than for giant-galax y-sized halos (the bottom 
panel). The result of iMoore et"all (|l999al ) (thick 
dashed curve) is similar to those for our rel atively 
subha lo-rich halos, while that of Dieman d et al.l 
( 20081 ) (thick dotted curve) is comparable to 
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Fig. 2. — Cumulative numbers of subhalos as a 
function of their maximum rotation velocities Vc 
normalized by those of the parent halos Vp. (a) 
68 galaxy-sized halos with 1.5 x IO^^Mq < M < 
3 X 1O^^M0 (top), (b) 57 giant-galaxy-sized ha- 
los with 3 X W'^Mq < M < 1 X lO^^^o (bot- 
tom). Three thick solid curves show the aver- 
age (middle) and ±la values (top and bottom). 
Thick dott ed and dashed curves are the fitting for- 
mula fromlDiemand et all ( 20081 ) and the result of 
Moore et al.l ( 1999a[ ) for a galaxy-sized halo, re- 
spectively. The thin dashed curve with open tri- 
angles denotes the nu mber of dwarf galaxies in our 
galaxy (|Matedll998l ). The open circles with error 
bars show t he number of dwarf g alaxies in the Lo- 
cal Group (jP'Onghia et al.ll2007l ). 



Fig. 3. — Cumulative function F of parent halos 
as a function of the number of subhalos with the 
normalized circular velocity larger than 0.1, A'^>o.i- 
The thick solid curve shows the result for galaxy- 
sized halos. The thin solid curve shows the re- 
sult for giant galaxy-sized halos. The thin dashed 
curve is for the result of Paper I. The arrows in the 
botto m show the results of s everal previous papers 
(D07 : iDiema nd et al.l l2q07l : p 08: 'Diemand et aL 



20081 K07: iKase et al 

mM). 



2007; M99: .Moore et al 



our average value. The thin dashed curve with 
open tri angles denote s the dwarf galaxies in our 
Galaxy (jMated 119981 ). and open circles with er- 
ror bars sh ow the dwarf galaxies in the Lo- 
cal Group (jP'Onghia et al.l l2007i l. Note that 
these observational results do not include satel- 



lites recently found by SPSS (IWill man et al| 
2005HBelokurov et al-lbOOllZucker et al. 2006a.b[ 



Belokurov et al.ll2007:llrwin et al.ll2007: Walsh et al 



20071 ) ■ and thus clearly are underestimates at 
least for the rotation velocity less than lOkms"^ 
{VjVp < 0.05). 

Figure [3] shows the cumulative distribution 
of the subhalo abundance of the parent halos. 
We used the number of subhalos with the nor- 
malized circular velocity larger than 0.1, A'^>o.i, 
as the measure of the subhalo abundance. We 
can see that the variation of the subhalo abun- 
dance of galaxy-sized halo is larger than that 
for giant-galaxy-sized halos. The abundances of 
the richest halos are roughly the same for both 
galaxy-sized and giant-galaxy-sized halos, but the 
distribution of galaxy-sized halos is extended to 
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Fig. 4. — Dependence of the subhalo abundance iV>o.io on (a) mass M12 in unit of (left), (b) 

concentration parameter c ^ ro/R^ (right). The white circles denote the halos with 1.5 x 10^'^ Mq < M < 
3 X lO^^MfT). The black circles denote the halos with 3 x lO^^Mm < Af < 1 x lO^^Mm. 



smaller values of A^>o.i- We also show the re- 
sult of sev e ral pr eviou s papers. The result s of 
Kase et all (|2007t ) and IPiemand et"al] (|2008l ) lie 
in the middle of the distribution of our halos, while 



the result of iMoore et al.l (|l999ar ) co rresponds to 
the richest ones. The differences of Moore et al 



1999aD . iKase et al.l (|2007^ . and iDiemand et al 



(|2008[ ) are probably p artly due to the d iffer- 



ence in the cosmo logy. iMoore et al. I (ll999a) and 



Kase et al.l (l2007l) adopted s t andar d CDM with- 
out A, while iDiemand et all (I2OO8I) used A CDM 
Our present calculation and Diema nd et al. I (|2008l) 
used quite similar ACDM cosmological param- 
eters. We will return to the question why the 
result of lDiemand et al.l (j2007l ) corresponds to our 
most subhalo-poor halos in Section 4.1. 

The dashed curve shows the result of Paper I, 
for galaxy-sized halos. The halos Paper I have 
somewhat smaller number of subhalos, but since 
there are only 10 galaxy-sized halos, the difference 
between Paper I result and the present result is 
well within the statistical variation (not significant 
for the KS test with 5% criterion). 

3.3. The cause of the variation of subhalo 
abundance 

Figure HKa) shows the dependence of the par- 
ent halo on the mass. As already shown in Fig- 
ure [21 the variation of the subhalo abundance of 
galaxy-sized halos is significantly larger than that 
for giant-galaxy-sized halos. 



Figure [3Kb) shows the dependence on the con- 
centration parameter c. Here, c = ro/i?v, where 
rn is the scale radius of the fit to the Moore profile 



(|Moore et al.lll999d ). and i?v is the virial radius 



of the halo. The number of subhalos shows fairly 
tight correlation with the concentration parame- 
ter. More centrally concentrated halos have fewer 
subhalos. This correlation might mean that the 
subhalo abundance depends on the halo formation 
epoch. Halos with large values of c are generally 
formed earlier. This difference in the formation 
time might be the reason for the difference in the 
subhalo abundance. 

In Paper I, we found dependences on the local 
density within 5Mpc sphere normalized by the av- 
erage density and also on the distance to the near- 
est influential larger-sized halo, which has three 
times or more mass of the parent halos. However, 
with the improved statistics of the present paper, 
we could not confirm these results. 

Figure [HI and [Sj show three subhalo-poor ha- 
los (top) and three subhalo-rich halos (bottom) 
at z = 6.04, z = 3.21, z = 1.04, and z = (left 
to right). They correspond to 0.85, 2.0, 5.68, and 
13.7Gyr. In Figure [6l only particles which lie in 
the halo virial radius at z = are plotted. As 
can be seen from the panels for z = in Fig- 
ure [H there is a large difference of the subhalo 
abundance. The large difference in the subhalo 
abundance is already clear at z = 1.04. At this 
moment, subhalo-poor halos have almost finished 
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Fig. 5. — Snapshots of three subhalo-poor halos (top) and three subhalo-rich halos (bottom). Starting from 
the left, epoch are z = 6.04, z = 3.21, z = 1.04, and z = 0, width are O.SMpc, 1.6Mpc, 1.6Mpc, and O.SMpc. 
The center of snapshot is defined as the center of mass of particles which lie in the halo virial radius at 
z = 0.0. 7 
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Fig. 6. — Same as Figure [5l However, only particles which lie in the halo virial radius at z = are plotted. 
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the assembling (except for the one in the top row, 
which still contains two large subhalos). On the 
other hand, subhalo-rich halos contain many more 
large halos which are in the process of first fall- 
back. Moreover, they are part of much larger, 
Mpc-scale filamentary structures. In panels for 
z = 3.21, again, the difference between subhalo- 
poor halos and rich halos is rather clear. Rich ones 
contain many small structures scattered all over 
the distribution of particles which form the halo 
at 2: = 0. In the case of poor halos, at z = 3.21 sin- 
gle high-density structure at the center is clearly 
visible. 

When we compare the snapshot at high-z in 
Figure we can also see that the physical size of 
the distribution of particles is significantly smaller 
for poor halos, compared to rich halos. Also, the 
shape of the halo is smoother and rounder for poor 
halos. 

Thus, from the way these halos are formed, as 
shown in figures [5] and [6l it is rather clear why 
there is large variation in the subhalo abundance. 
Poor halos are formed from more centrally concen- 
trated initial condition and formed earlier, while 
rich halos are formed from less concentrated struc- 
ture and formed later. 

In order to quantify the difference in the way 
the halos are formed, we tried many different mea- 
sures, such as the time of the last major merger 
and the time at which the mass within the virial 
radius (or its certain fraction) exceeds 50% (or 
other fraction) of the final mass. Most measures 
we tried show the correlation with the subhalo 
abundance, but that was weaker than that of the 
concentration parameter c. One measure with a 
fair success is the half-mass radius at the maxi- 
mum expansion, i?h,max, defined as the maximum 
value of the radius in which the half of the final 
virial mass is enclosed. To calculate this radius, 
we used the potential minimum as the center of 
the sphere. 

Figure [7| shows the evolution of the half-mass 
radius The solid and dashed curves show the 
results of the rich and poor halos, respectively. 

We can see that i?h of poor halos reach the 
maximum values much earlier than those for rich 
halos, and the maximum values are smaller for 
poor halos. Also, _Rh of poor halos seem to settle 
down faster than those of subhalo-rich halos. 



Figure [8] shows the relation between the concen- 
tration parameter c and the maximum half-mass 
radius -Rh,max normalized by the halo virial radius 
i?v (left), and also A^>o.i and i?h,max- We can see 
that the correlation is fairly tight for both cases. 
Thus, the number of subhalos is at least partly 
determined by the way the halo is assembled. 

Figure [QlfTTl show three subhalo-poor halos 
(top) and three subhalo-rich halos (bottom) at 
z = 6.04, z = 3.21, z = 1.04, and z = (left 
to right). Their i?h,max/^v are almost same (up- 
per to bottom, 0.71, 0.74, 0.68, 0.69, 0.78, 0.81) 
and their concentration parameter are also rather 
close (7.0, 7.0, 4.8, 6.1, 6.1, 5.9). Two of them are 
halos with 1.5 x W^Mq < M < 2.0 x lO^^Mg 
(third and sixth from top). Others are halos with 
2.0 X 10^2^/0 <M< 2.5 X 10^2^0. However, we 
can see the large difference in the subhalo abun- 
dance from the panels for z = in Figure 1111 
As can be seen from the panels for z — 1.04, the 
difference of the shape of the halo is clear. At 
this moment, subhalo-poor halos have one or two 
large subhalos. By contrast, subhalo-rich halos 
have many small subhalos. In panels for high-z, 
we can see the same feature for subhalo-poor and 
subhalo-rich halos as seen in Figure [HI 

These mean that subhalo-poor halos with large 
maximum half-mass radius are formed from more 
centrally concentrated initial condition, formed 
earlier and have a small number of large compan- 
ions. In such a halo, the maximum half-mass ra- 
dius is large because of these companions. On the 
other hand, subhalo-rich halos have many small 
subhalos at z = 1.04 and many merger events af- 
ter that. This might be the reason for the large 
difference in the subhalo abundance of halos with 
same maximum half-mass radius. 

4. Discussions and Summary 

4.1. Comparison with Aquarius simula- 
tion 



Springel et al. ( 2008h claimed that they found 
very small halo-to-halo variations in the subhalo 
abundance. They used the definition of the sub- 
halo count different from what was used in our 
work and most previous works. In this section, we 
compare our result with theirs. 

Figure[T2]shows the cumulative numbers of sub- 
halos as a function of their maximum rotation ve- 
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Fig. 7. — Evolution of the halo half- mass radius Rh{t). The solid curves show the results of three subhalo- 
rich halos. The dashed curves show the results of three subhalo-poor halos. The left panel is for halos with 
1.5 X lO^^M© < M < 2.0 X lO^^M©. The right panel is for halos with 2.0 x 10^^^^ < M < 2.5 x 10^^^^ 
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Fig. 8. (a) Concentration c plottc;(l against the maximum half-mass radius i?h.max normalized by the 
halo virial radius Ry (left), (b) The subhalo abundance A^>o.io plotted against the maximum half-mass 
radius i?h,max normalized by the halo virial radius Ry (right). Black circles, white circles, and crosses 
denote the halos with 2.5 x < Af < 3 x 10^^^^^ 2.O x 10^^^^ < M < 2.5 x W^^Mq, and 

1.5 X IO^^Mq < M < 2.0 X respectively. 
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Fig. 9. — Snapshots of three subhalo-poor halos (top) and three subhalo-rich halos (bottom). Starting from 
the left, epoch are z = 6.04, z = 3.21, z = 1.04, and z = 0, width are O.SMpc, 1.6Mpc, 1.6Mpc, and O.SMpc. 
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Fig. 11. — Same as Figure[51 However, only particles which lie in the halo virial radius at z = are plotted. 
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locities Vc normalized by the circular velocity of 
the parent halos at r^o- We define the radius rso 
as inside which the averaged spherical overdensity 
is 50 times the critical value and count the num- 
ber of subhalos within rso. This d efinition is the 
same as that used for Figure 10 in lSpringel et al. 
(|2008l ). We plot the result of the 44 galaxy-sized 
halos with M50 < 3.0 x lO^^M©, where M50 is the 
mass within r^Q. 

If we measure the subhalo abundance again us- 
ing iV>o.i, the number of subhalos with normal- 
ized rotation velocity m ore than 0.1, w e can c on- 
clude that our result and Springel et al. ( 2008h re- 
sult are practically identical. However, our re- 
sults for iV>o.i may be affected by the resolu- 
tion, because V50 is typically 20% smaller than 
the maximum rotation velocity Vp. Thus, our ha- 
los are p robably slightl y more subhalo-rich than 
those of Springel et al.l ( 2008h . but the halo-to- 
halo variation , calcu lated using the definition of 
Springel et al.l (|2008r ). is largely similar, though 
it is difficult to tell much on the variation in 



Springel et al.l (j2008l ) result because of small num- 
ber of runs. We show three results here, and the 
total number of runs performed was six. 

Compared with Figure [2l the variation of the 
subhalo abundance in Figure [12] is significantly 
smaller. This result is quite natural because of 
the two differences in the definition of the subhalo 
count. First, the use of r^Q means that many sub- 
halos are at the outermost region which are not 
virialized. Though subhalos in this outer region 
are gravitationally bound to the main halo, their 
evolution is not affected by the structure of the 
main halo simply because they are far from the 
central high-density region. Second, the use of 
V50 instead of Vp would also reduce the halo-to- 
halo variation, since centrally concentrated halos 
have higher Vp. 

Springel et al. ( 2008f l claimed that the small 
halo-to-halo variation they found contradicted 
with the result of Paper I. Their rms halo-to- 
halo scatter was about 8%, while ours (in Figure 
2) is about 41%. From the comparison above, it 
is clear that this "contradiction" they found is 
partly due to the difference in the definition of the 
subhalo abundance. The halo-to-halo scatter as 
seen in Figure [12] is smaller than one in Figure [2] 
However, there is still some difference between our 
results and their ones. The remaining difference 
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Fig. 12. — Cumulative numbers of subhalos as a 
function of their maximum rotation velocities Vc 
normalized by the circular velocity of the parent 
halos at r^Q. Three thick solid curves show the 
average (middle) and ±la values (top and bot- 
tom) . The thick dash ed curves show the results of 
ISpringel etaP l|2008L simulation Aq-A-1, Aq-C-2, 
Aq-E-2). 
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might be caused by their way to select parent ha- 
los. They selected the halos which host normal 
spirals in their semianalytic modelling. 



Th e systematic difference between lSpringel et al, 
(j2008) results and our re sult may be due to the 



wav lSpringel et al.l (|2008f ) selected their halos for 
high-resolution calculations. They selected halos 
without massive neighbors at z — 0, and they also 
put a criterion that in the semi-analytic calcula- 
tion the halo should host a late-type galaxy. This 
would mean the selected halos should not experi- 
ence major merger events at low-z, and probably 
implies that their formation epoch is early. Thus, 
it is natural that their halos have systematically 
smaller numbers of subhalos than ours have. They 
found systematic difference between their re sult 
and Via Lactea (jPiemand et al.l l2007i |2008|) re- 



sults. However, our result agrees with both. 

4.2. Is the missing-dwarf problem solved? 

The difference with the observed number of 
dwarf galaxies in Local Group or our galaxy is 
around a factor of 2 for the most subhalo-poor ones 
in our simulation, at the rotation velocity 10% of 
that of the parent halo. For smaller velocities, the 
difference is larger, but recent SDSS results indi- 
cated that there are still many dwarf galaxies not 
identified. Thus, the missing dwarf problem does 
not seem to be a very serious problem. A factor- 
of-two difference may be easily explained by the ef- 
fect of baryons. Also, Springel et al. ( 20081 ) found 
the subhalo abundance even lower than the lowest 
number we found. Thus, if our galaxy corresponds 
to most subhalo-poor halos formed in ACDM uni- 
verse, there is little discrepancy between the ob- 
served number of dwarf galaxies and the expected 
number of dark matter subhalos. 

Obvious question is why our galaxy belongs 
to most subhalo-poor halos, but as we found the 
halos with high central concentration and early 
formation epoch have generally smaller number 
of subhalos. Also, since our galaxy is a spiral 
galaxy, it is unlikely that it experienced recent ma- 
jor merger event. If we look at the Local Group as 
a whole, it will experience a major merger in some 
future, but in this case we should look at two main 
subhalos, our galaxy and M31, since their struc- 
ture determine the total number of subhalos in the 
Local Group. 



4.3. Summary 

We simulated 1600'^ particles with mass of 1.0 x 
1O^M0 in a 46.48Mpc cubic box to create an unbi- 
ased sample of halos with resolution high enough 
to determine the abundance of subhalos. We con- 
firmed the large variation of subhalo abundance 
we found in Paper I, but could not confirm the 
dependence on environmental parameters such as 
average density or distance to the nearest mas- 
sive halo. However, we found clear dependence 
on the concentration parameter, which means the 
number of subhalos is determined primarily by the 
initial condition and formation time. 

The "mis s ing-dw arf problem" pointed out by 
Moore et al.l (|l999al ) turned out to be much less 
serious, since now the discrepancy is around a fac- 
tor of 2. A detailed study of baryon evolution of 
dwarf galaxies might be necessary to close up this 
remaining discrepancy. 
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